home
***
CD-ROM
|
disk
|
FTP
|
other
***
search
/
The PC-SIG Library 10
/
The PC-Sig Library - Shareware for the IBM PC and Compatibles (PC-SIG)(Tenth Edition Disks 1-2804)(1991).iso
/
PC_SIGCD
/
09
/
2
/
DISK0921.ZIP
/
TNODE.BAS
< prev
next >
Wrap
BASIC Source File
|
1987-09-29
|
12KB
|
454 lines
' PROGRAM "TNODE"
' COPYRIGHT (C) 1983 BY DAVID EAGLE
' IBM-PC << QUICKBASIC COMPILER VERSION 3.0 >>
' PUBLIC DOMAIN FOR IBM-PC ON NOVEMBER 20, 1986
' COMPUTES ASCENDING NODE CROSSINGS OF EARTH SATELLITES
' DATE OF CROSSING; MONTH, DAY, YEAR
' GMT OF CROSSING; HOURS, MINUTES, SECONDS
' WEST LONGITUDE OF CROSSING; DEGREES
'*********************************************************
DEFDBL A-Z
DIM SHARED MONTH$(12)
COMMON SHARED ECC,INC,ARGPER,TN,THOUR,TMIN,TSEC,CDATE$
COMMON SHARED MONTH%,DAY%,YEAR%,JDATE,WLONG,CTIME$
CONST PI=3.141592653589793D0
CONST PI2=2D0*PI
CONST PIDIV2=.5D0*PI
CONST DTR=PI/180D0
CONST RTD=180D0/PI
CONST MU=1.407645794D16
CONST REQ=20925656.2D0
CONST J2=.0010828D0
CONST WE=7.29211515D-5
MONTH$(1)="January"
MONTH$(2)="February"
MONTH$(3)="March"
MONTH$(4)="April"
MONTH$(5)="May"
MONTH$(6)="June"
MONTH$(7)="July"
MONTH$(8)="August"
MONTH$(9)="September"
MONTH$(10)="October"
MONTH$(11)="November"
MONTH$(12)="December"
'**********************************************************
CLS
PRINT
PRINT
PRINT "Program TNODE"
PRINT "(C) Copyright 1982 by David Eagle
PRINT
PRINT "Microsoft QuickBASIC Compiler"
PRINT "(C) Copyright Microsoft Corp. 1982-1987"
PRINT
CALL KEYCHECK
CLS
PRINT
PRINT
INPUT "Introduction ( y = yes, n = no )";A$
IF INSTR("yY",A$) THEN CALL INTRO
' PROMPT USER FOR INPUTS
CLS
PRINT
PRINT
PRINT TAB(12);"Program TNODE"
PRINT
SATELLITE:
PRINT
PRINT "Orbital period ( minutes )"
INPUT TKEPLER
TKEPLER=TKEPLER*60
PRINT
PRINT "Orbital inclination ( degrees [ 0 - 180 ] )"
INPUT INC
INC=INC*DTR
PRINT
PRINT "Orbital eccentricity ( non-dimensional )"
INPUT ECC
IF ECC>0D0 THEN
PRINT
PRINT "Argument of perigee ( degrees [ 0 - 360 ] )"
INPUT ARGPER
ARGPER=ARGPER*DTR
ELSE
ARGPER=0D0
END IF
PRINT
REVENT:
CLS
PRINT
PRINT
PRINT TAB(16);"Reference Event"
PRINT
PRINT
PRINT TAB(10);"< 1 > Equatorial crossing"
PRINT
PRINT TAB(10);"< 2 > Apogee"
PRINT
PRINT
PRINT "Selection"
INPUT TYPE%
PRINT
PRINT
PRINT "Date of reference event ( month [ 1 - 12 ], day [ 1 - 31 ], year [ YYYY ] )"
PRINT "< For example, October 21, 1985 is input as 10,21,1985 >"
INPUT MONTH%,DAY%,YEAR%
PRINT
PRINT "GMT of reference event"
PRINT "( hours [ 0 - 23 ], minutes [ 0 - 59 ], seconds [ 0 - 59 ] )"
PRINT "< For example, 8:45:30 p.m. is input as 20,45,30 >"
INPUT THOUR,TMIN,TSEC
PRINT
PRINT "West longitude of reference event ( degrees [ 0 - 360 ] )"
INPUT WLONG
WLONG=WLONG*DTR
PRINT
PRINT "Number of crossings to compute"
INPUT NCROSS%
CLS
PRINT
PRINT TAB(8);"CALENDAR";TAB(31);"GMT";TAB(50);"WEST";
PRINT TAB(10);"DATE";TAB(31);"TIME";TAB(48);"LONGITUDE"
PRINT
' COMPUTE REFERENCE JULIAN DATE AND TIME BETWEEN CROSSINGS
JDATE=367D0*YEAR%-INT(7*((YEAR%+INT((MONTH%+9)/12))/4))+INT(275*MONTH%/9)+DAY%+1721014D0
A=TKEPLER/PI2
A=((MU*A*A)^(1D0/3D0))/REQ
B=1D0+ECC*COS(ARGPER)
C=1D0-ECC*ECC
D=A*A
E=B*B
F=SIN(INC)
G=COS(INC)
H=1D0/(A*C)
TN=TKEPLER*(1D0-.5D0*J2*((4D0-5D0*F*F)/(D*SQR(C)*E)+2D0*E*B/(C*C*C*D)))
RDOT=-1.5D0*J2*PI2*H*H*G/TKEPLER
DLONG=TN*(WE-RDOT)
A2=TN/3600D0
DHOUR=INT(A2)
A4=60D0*(A2-DHOUR)
DMIN=INT(A4)
DSEC=INT(60D0*(A4-DMIN)+.5D0)
IF TYPE%=2 THEN CALL APOGEE
WLONG=WLONG-DLONG
THOUR=THOUR-DHOUR
TMIN=TMIN-DMIN
TSEC=TSEC-DSEC
' COMPUTE SUCCESSIVE CROSSINGS
FOR I%=1 TO NCROSS%+1
WLONG=WLONG+DLONG
IF WLONG>PI2 THEN WLONG=WLONG-PI2
THOUR=THOUR+DHOUR
TMIN=TMIN+DMIN
TSEC=TSEC+DSEC
CALL CTIME
CALL GDATE
CALL PRTSUB
NEXT I%
PRINT
CALL KEYCHECK
CLS
PRINT
PRINT
INPUT "Another selection ( y = yes, n = no ) ";A$
IF INSTR("nN",A$) THEN END
PRINT
INPUT "Another satellite ( y = yes, n = no ) ";A$
IF INSTR("yY",A$) THEN GOTO SATELLITE
PRINT
INPUT "Another reference event ( y = yes, n = no ) ";A$
IF INSTR("yY",A$) THEN GOTO REVENT
END
'**********************************************************
SUB GDATE STATIC
' GREGORIAN DATE SUBROUTINE
A=INT((JDATE-1867216.25D0)/36524.25D0)
A=JDATE+A-INT(A/4D0)+1D0
B=A+1524D0
C=INT((B-122.1D0)/365.25D0)
D=INT(365.25D0*C)
E=INT((B-D)/30.6001D0)
DAY%=B-D-INT(30.6001D0*E)
IF E<13.5 THEN
MONTH%=E-1
ELSE
MONTH%=E-13
END IF
IF MONTH%>2.5 THEN
YEAR%=C-4716
ELSE
YEAR%=C-4715
END IF
CDATE$=MONTH$(MONTH%)+STR$(DAY%)+","+STR$(YEAR%)
END SUB
'**********************************************************
SUB PRTSUB STATIC
' PRINT SUBROUTINE
PRINT TAB(5);CDATE$;
PRINT ;TAB(26);CTIME$;
PRINT USING "###.##";TAB(48);WLONG*RTD
END SUB
'**********************************************************
SUB KEYCHECK STATIC
' CHECK USER RESPONSE SUBROUTINE
PRINT
PRINT TAB(15);"< press any key to continue >";
A$=""
WHILE A$=""
A$=INKEY$
WEND
END SUB
'**********************************************************
SUB INTRO STATIC
' INTRODUCTION SUBROUTINE
CLS
PRINT TAB(10);" INTRODUCTION "
PRINT
PRINT TAB(10);" PROGRAM 'TNODE' IS AN INTERACTIVE "
PRINT TAB(10);"BASIC COMPUTER PROGRAM WHICH CAN BE USED"
PRINT TAB(10);"TO DETERMINE ASCENDING NODE CROSSINGS OF"
PRINT TAB(10);"EARTH SATELLITES. THE ASCENDING NODE IS "
PRINT TAB(10);"THE INTERSECTION OF THE EARTH'S EQUATORIAL"
PRINT TAB(10);"PLANE AND THE SATELLITE ORBIT PLANE AS"
PRINT TAB(10);"THE SATELLITE 'ASCENDS' FROM SOUTH TO"
PRINT TAB(10);"NORTH. THE TIME AND GEOGRAPHIC LONGITUDE"
PRINT TAB(10);"OF AN ASCENDING NODE CAN BE USED TO "
PRINT TAB(10);"PREDICT THE SATELLITE'S POSITION AT OTHER"
PRINT TAB(10);"TIMES RELATIVE TO AN OBSERVER ANYWHERE"
PRINT TAB(10);"IN THE WORLD. "
PRINT
PRINT TAB(10);" USER INPUTS AND SELECTIONS "
PRINT
PRINT TAB(10);" PROGRAM 'TNODE' WILL PROMPT THE USER"
PRINT TAB(10);"FOR SEVERAL INPUTS NECESSARY FOR THE "
PRINT TAB(10);"SOFTWARE TO WORK PROPERLY. THIS IS A "
PRINT TAB(10);"DESCRIPTION OF THESE REQUESTS AND A "
PRINT TAB(10);"DISCUSSION OF THE PROPER USER RESPONSE. "
CALL KEYCHECK
CLS
PRINT
PRINT TAB(10);"Orbital period ( minutes )"
PRINT
PRINT TAB(10);"THE USER SHOULD RESPOND WITH THE ORBITAL"
PRINT TAB(10);"PERIOD OF THE SATELLITE IN MINUTES. "
PRINT
PRINT TAB(10);"Orbital inclination ( degrees )"
PRINT
PRINT TAB(10);"THE RESPONSE TO THIS REQUEST SHOULD BE "
PRINT TAB(10);"THE ORBITAL INCLINATION OF THE SATELLITE"
PRINT TAB(10);"IN DEGREES. THIS INPUT SHOULD BE A "
PRINT TAB(10);"NUMBER BETWEEN 0 AND 180 DEGREES. "
PRINT
PRINT TAB(10);"Orbital eccentricity ( non-dimensional )"
PRINT
PRINT TAB(10);"THE USER SHOULD INPUT THE ECCENTRICITY "
PRINT TAB(10);"OF THE ORBIT. THIS IS A NON-DIMENSIONAL "
PRINT TAB(10);"NUMBER BETWEEN 0 AND 1. "
PRINT
PRINT TAB(10);"Argument of perigee ( degrees )"
PRINT
PRINT TAB(10);"THE USER SHOULD RESPOND WITH THE ORBIT'S"
PRINT TAB(10);"ARGUMENT OF PERIGEE IN DEGREES. "
CALL KEYCHECK
CLS
PRINT
PRINT TAB(10);" Reference Event "
PRINT
PRINT TAB(10);" <1> Equatorial crossing "
PRINT
PRINT TAB(10);" <2> Apogee "
PRINT
PRINT TAB(10);"Selection? "
PRINT
PRINT TAB(10);"HERE THE USER SHOULD INDICATE THE TYPE "
PRINT TAB(10);"OF REFERENCE EVENT. MOST PREDICTIONS ARE"
PRINT TAB(10);"BASED ON EQUATORIAL CROSSINGS WHILE SOME"
PRINT TAB(10);"PREDICTIONS, SUCH AS THOSE FOR OSCAR 10,"
PRINT TAB(10);"ARE BASED ON APOGEE. "
PRINT
PRINT TAB(10);"Date of reference event ( month, day, year )?"
PRINT
PRINT TAB(10);"AT THIS POINT THE USER SHOULD INPUT THE "
PRINT TAB(10);"NUMERICAL EQUIVALENT OF THE DATE OF THE "
PRINT TAB(10);"REFERENCE EVENT. FOR EXAMPLE, MAY 12, "
PRINT TAB(10);"1982 IS INPUT AS '5,12,1982'. "
CALL KEYCHECK
CLS
PRINT
PRINT TAB(10);"GMT of reference event ( hours, minutes, seconds )?"
PRINT
PRINT TAB(10);"THE USER SHOULD RESPOND WITH THE GREENWICH"
PRINT TAB(10);"MEAN TIME (GMT), IN HOURS, MINUTES AND "
PRINT TAB(10);"SECONDS CORRESPONDING TO THE REFERENCE "
PRINT TAB(10);"EVENT DESCRIBED ABOVE. "
PRINT
PRINT TAB(10);"West longitude of reference event ( degrees )?"
PRINT
PRINT TAB(10);"THE RESPONSE TO THIS REQUEST SHOULD BE "
PRINT TAB(10);"THE WEST LONGITUDE, IN DEGREES, OF THE "
PRINT TAB(10);"REFERENCE EVENT. WEST LONGITUDE IS EQUAL"
PRINT TAB(10);"TO 360 DEGREES MINUS EAST LONGITUDE. "
PRINT
PRINT TAB(10);"Number of crossings to compute? "
PRINT
PRINT TAB(10);"THE USER SHOULD RESPOND WITH THE NUMBER "
PRINT TAB(10);"OF ADDITIONAL CROSSINGS THE PROGRAM "
PRINT TAB(10);"SHOULD COMPUTE. DATA ABOUT THESE CROSSINGS"
PRINT TAB(10);"WILL BE BASED ON THE REFERENCE EVENT"
PRINT TAB(10);"DESCRIBED ABOVE. "
CALL KEYCHECK
CLS
PRINT
PRINT TAB(10);" AFTER THE PROGRAM HAS RUN, IT WILL ASK"
PRINT TAB(10);"THE USER FOR ANOTHER SELECTION. THE USER"
PRINT TAB(10);"SHOULD RESPOND WITH 'Y' IF HE OR SHE "
PRINT TAB(10);"DESIRES THE SELECTION OR 'N' IF NOT. "
PRINT
PRINT TAB(10);"Another selection ( y = yes, n = no )? "
PRINT
PRINT TAB(10);"THE USER SHOULD RESPOND WITH 'N' TO EXIT"
PRINT TAB(10);"THE PROGRAM. "
PRINT
PRINT TAB(10);"Another satellite ( y = yes, n = no )? "
PRINT
PRINT TAB(10);"THE USER SHOULD SELECT 'Y' TO COMPUTE "
PRINT TAB(10);"CROSSINGS FOR ANOTHER SATELLITE. "
PRINT
PRINT TAB(10);"Another reference event ( y = yes, n = no )? "
PRINT
PRINT TAB(10);"THE USER SHOULD INPUT 'Y' TO COMPUTE "
PRINT TAB(10);"EQUATORIAL CROSSINGS BASED ON ANOTHER "
PRINT TAB(10);"REFERENCE EVENT. "
CALL KEYCHECK
CLS
PRINT
PRINT TAB(10);" PROGRAM OUTPUT"
PRINT
PRINT TAB(10);"PROGRAM 'TNODE' WILL OUTPUT THE CALENDAR"
PRINT TAB(10);"DATE, GMT AND WEST LONGITUDE OF ASCENDING"
PRINT TAB(10);"NODE CROSSINGS. GMT IS PRINTED IN HOURS,"
PRINT TAB(10);"MINUTES AND SECONDS AND THE WEST LONGITUDE"
PRINT TAB(10);"IS PRINTED IN DECIMAL DEGREES."
CALL KEYCHECK
END SUB
'**********************************************************
SUB APOGEE STATIC
' APOGEE EVENT SUBROUTINE
D=PI2-ARGPER
A=SIN(D)*SQR(1D0-ECC*ECC)
B=ECC+COS(D)
CALL ATAN3(A,B,C)
D=PI-C+ECC*SIN(C)
TREL=-.5D0*D*TN/PI
IF TREL<0D0 THEN TREL=TN+TREL
A=TREL/3600D0
THOUR=THOUR+INT(A)
B=(A-INT(A))*60D0
TMIN=TMIN+INT(B)
TSEC=TSEC+INT((B-INT(B))*60D0+.5D0)
CALL CTIME
ARGLAT=ARGPER+PI
A=SIN(ARGLAT)*COS(INC)
B=COS(ARGLAT)
CALL ATAN3(A,B,C)
WLONG=WLONG+C+(WE+RDOT)*TREL
IF WLONG>=PI2 THEN WLONG=WLONG-PI2
END SUB
'**********************************************************
SUB ATAN3(A,B,C) STATIC
' FOUR QUADRANT ARC TANGENT SUBROUTINE
IF ABS(A)<1D-8 THEN
C=(1D0-SGN(B))*PIDIV2
EXIT SUB
ELSE
C=(2D0-SGN(A))*PIDIV2
END IF
IF ABS(B)<1D-8 THEN
EXIT SUB
ELSE
C=C+SGN(A)*SGN(B)*(ABS(ATN(A/B))-PIDIV2)
END IF
END SUB
'**********************************************************
SUB CTIME STATIC
' TIME SUBROUTINE
IF TSEC>=60D0 THEN
TSEC=TSEC-60D0
TMIN=TMIN+1D0
END IF
IF TMIN>=60D0 THEN
TMIN=TMIN-60D0
THOUR=THOUR+1D0
END IF
IF THOUR>=24D0 AND TMIN>0D0 THEN
THOUR=THOUR-24D0
JDATE=JDATE+1D0
PRINT
END IF
CTIME$=STR$(THOUR)+" h"+STR$(TMIN)+" m"+STR$(TSEC)+" s"
END SUB